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Temperature inversions occur in nature, e.g., in the solar corona and in interstellar molecular 
clouds: somewhat counterintuitively, denser parts of the system are colder than dilute ones. We 
propose a simple and appealing way to spontaneously generate temperature inversions in systems 
with long-range interactions, by preparing them in inhomogeneous thermal equilibrium states and 
then applying an impulsive perturbation. In similar situations, short-range systems would typically 
relax to another thermal equilibrium, with uniform temperature profile. By contrast, in long-range 
systems, the interplay between wave-particle interaction and spatial inhomogeneity drives the system 
to nonequilibrium stationary states that generically exhibit temperature inversion. We demonstrate 
this mechanism in a simple mean-field model and in a two-dimensional self-gravitating system. Our 
work underlines the crucial role the range of interparticle interaction plays in determining the nature 
of steady states out of thermal equilibrium. 
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Stationary states far from thermal equilibrium occur 
in nature. In some cases, e.g., in the solar corona and 
in interstellar molecular clouds, such states exhibit tem¬ 
perature inversion: denser parts of the system are colder 
than dilute ones. This work is motivated by an attempt 
to explain how such a counterintuitive effect may sponta¬ 
neously arise in nonequilibrium states, unveiling its min¬ 
imal ingredients and the underlying physical mechanism. 
To this end, we start with asking a simple yet physically 
relevant question: What happens if an isolated macro¬ 
scopic system in thermal equilibrium is momentarily dis¬ 
turbed, e.g., by an impulsive force or a “kick”? If the 
interactions among the system constituents are short- 
ranged, collisions redistribute the kick-injected energy 
among the particles, yielding a fast relaxation to a new 
equilibrium, with a Maxwellian velocity distribution and 
a uniform temperature across the system. 

Is the scenario the same if instead the interactions are 
long-ranged [1]? For long-range systems, collisional ef¬ 
fects act over a characteristic time Tcoii that, unlike short- 
range systems, scales with the system size N, diverging as 
—>■ oo [2]. Asa result, a macroscopic system with long- 
range interactions starting from generic initial conditions 
will attain thermal equilibrium only after extremely long 
times, often exceeding typical observation times. Exam¬ 
ples of long-range systems are self-gravitating systems, 
for which, e.g., Tcoh — 10^° years for globular clusters and 
orders of magnitude larger than the age of the universe 
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for galaxies [3, 4]. The collisionless evolution of long- 
range interacting systems for times shorter than Tcou is 
governed by the Vlasov (or collisionless Boltzmann) equa¬ 
tion [2] . When kicked out of thermal equilibrium, a long- 
range interacting system relaxes to a Vlasov-stationary 
state, and thermal equilibrium is just one of infinitely 
many possible states. Predicting which Vlasov state fol¬ 
lows a given initial condition is an extremely difficult 
problem [5], first addressed in [6] and still unsolved, de¬ 
spite some remarkable achievements for particular sys¬ 
tems and special initial conditions (see [7, 8]). 

Let us then ask a simpler question: How different from 
an equilibrium state is the stationary state the system 
relaxes to after the kick? Are there ways to characterize 
it, e.g., by unveiling some of its general features? The 
answer is yes; in this paper, we argue that, provided it 
started from a spatially inhomogeneous equilibrium state 

[9] , the system after the kick relaxes to a state with a non- 
uniform temperature profile. In short-range systems, by 
contrast, a non-uniform temperature profile may only oc¬ 
cur when the system is actively maintained out of equilib¬ 
rium, e.g., by a boundary-imposed temperature gradient, 
to counteract collisional effects. Remarkably, in a long- 
range system, the relaxed state after the kick generically 
exhibits temperature inversion, as we will explicitly show. 

As recalled above, temperature inversions are observed 

[10] in interstellar molecular clouds [11-13] and especially 
in the solar corona, where temperatures around 10® K 
that are three orders of magnitude larger than the tem¬ 
perature of the photosphere are attained [14]. Despite 
recent advances [15], the mechanism of coronal heating 
is not completely understood and remains one of the 
most important open problems in astrophysics [16]. Most 
attempts to explain such a phenomenon involve mecha- 
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nisms that actively inject energy [17] into the less dense 
regions of the system. A different possibility, suggested 
by Scudder [18-20], is referred to as velocity filtration; 
see also Ref. [21] and Appendix A. Consider a system 
of particles acted upon by an external field whose po¬ 
tential energy increases with height above a base level. 
Only particles with a sufficiently large kinetic energy can 
climb the potential well and reach a given height. If 
the velocity distribution at the base is Maxwellian, it re¬ 
mains as such with the system maintaining the same tem¬ 
perature at all heights in the stationary state. Instead, 
if the distribution is suprathermal, i.e., with tails fat¬ 
ter than Maxwellian, the temperature in the stationary 
state increases with height, with concomitant decrease of 
the density. It is argued [16, 22, 23] that velocity filtra¬ 
tion might not be the (only) mechanism behind coronal 
heating, but nevertheless provides a simple and appeal¬ 
ing explanation of how a counterintuitive temperature 
inversion occurs without steady energy injection in less 
dense parts of the system. 

Scudder’s original model neglects interparticle interac¬ 
tions and requires an “active” ingredient, i.e., an out-of¬ 
equilibrium suprathermal velocity distribution imposed 
as a boundary condition. It was recently shown [21] that 
temperature inversion occurs also in strongly interacting 
systems, provided the interactions are long-ranged, when 
a velocity distribution with suprathermal tails is given 
just as the initial condition of the dynamical evolution. 
Although much weaker than a nonthermal boundary con¬ 
dition, the latter is still an ad hoc requirement. How¬ 
ever, it is not necessary at all: in this paper, we demon¬ 
strate that temperature inversion emerges spontaneously 
in the stationary state reached after a long-range system 
is brought out of equilibrium by a perturbation acting for 
a very short time, and there is no need for a suprathermal 
initial distribution. Our work thus suggests that temper¬ 
ature inversions observed in astrophysical systems may 
be examples of a more general phenomenon, whose roots 
are in the long-range nature of the interparticle interac¬ 
tions. 

Let us consider the very general setting of a system of 
N interacting particles of mass m in d dimensions with 
Hamiltonian 

N 2 N N 
i=l i=l j>i 

where pi = jpij are the momenta conjugated to the po¬ 
sitions Ti, the potential energy is long-ranged, V{r) oc 
r““;0 < Of < d as r —>■ oo. For times t <C Tcoii, the 
dynamics is described by the single-particle phase space 
distribution /(r,p,t) obeying the Vlasov equation 

df 

^ + p-Vr/-Vrw[/]-Vp/ = 0, (2) 

where 

■«[/] = J drdpf{r,p,t)V{r) (3) 


is the mean-field potential energy. An initial condition 
/o = /(r, p, 0) chosen to be a stationary solution of Eq. 
(2) does not evolve in time. Otherwise, if the initial dis¬ 
tribution is not a stationary solution of (2), after a short 
transient (often referred to as “violent relaxation”, af¬ 
ter Lynden-Bell [6]), the system settles into a stable sta¬ 
tionary solution of Eq. (2) called a quasi-stationary state 
(QSS), in which the system remains trapped until, at 
t ~ Tcoii, collisional effects neglected in Eq. (2) drive the 
system towards thermal equilibrium [2]. To construct a 
first representative example to demonstrate our claims, 
we assume periodic coordinates so that boundary effects 
may be neglected [24] and expand the interparticle po¬ 
tential (that by definition is even) in a cosine Eourier 
series. We then set d = 1 and retain just the first Eourier 
term. The resulting model, the so-called Hamiltonian 
Mean-Field (HMF) model, describes a system of glob¬ 
ally interacting point particles moving on a circle with 
Hamiltonian [25] 

N 2 . N N 

Hhmf = 

2=1 2 = 1 j>i 

where idi € (—tt, tt] is the angular coordinate of the iih 
particle on the circle, while pi is the conjugated momen¬ 
tum; we have further assumed to = 1 and the interaction 
to be attractive, fixing the energy scale to unity, and 
scaled it by 1/N to make it extensive (Kac prescription, 
[2]). The Hamiltonian (4) is invariant under the 0(2) 
symmetry group. In thermal equilibrium, for energy den¬ 
sity £ = E/N smaller than Sf. = 3/4, the symmetry is 
spontaneously broken to result in a clustered state. The 
order parameter of clustering is the average magnetiza¬ 
tion [26] (TOa,,TOy) = jf (SfciCOSi?i,X]fciSindiy The 
HMF model is a simple system which, besides serving as 
a framework to study statics and dynamics of long-range 
systems, actually models physical systems like gravita¬ 
tional sheet models and free-electron lasers [2]. 

In order to study what happens when we “kick” the 
HMF system out of equilibrium, we studied its dynami¬ 
cal evolution via molecular dynamics (MD) simulations, 
performed by integrating [27] the equations of motion de¬ 
rived from the Hamiltonian (4). We prepared the system 
in an equilibrium state with = toq and ruy = 0, and 
with a Maxwellian velocity distribution corresponding to 
the equilibrium temperature. We let the system evolve 
until t = tg > 0, and then kicked it out of equilibrium 
by applying during a short time r an external magnetic 
field h along the x direction; thus, for to < t < to + r, the 
Hamiltonian (4) is augmented by Hh = —cosi?i. 
Here we present results for to = 100, r = 1, h = 10, 
and Too = 0.521, corresponding to an initial equilibrium 
temperature T = 0.4244. We considered up to V = 10^ 
particles. After the kick the magnetization starts oscillat¬ 
ing, and after a transient the oscillations damp down and 
a new stationary value m* < toq of the magnetization is 
reached. A typical time evolution of the magnetization is 
shown in Fig. 1. The fact that to* < toq is not surprising 
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FIG. 1. (Color online) HMF: time evolution of the magneti¬ 
zation m with N = 10^ (solid red line). Inset: m{t) compared 
to the equilibrium value meq = 0 (dotted black line) for longer 
times. 
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FIG. 2. (Color online) HMF: temperature profile r(i9) (blue 
solid line) and density profile n(i9) (red dashed line), measured 
in the QSS at t = 10"*. 


because the system gains energy during the interaction 
with the external field h, resulting in an energy density 
e* > Ec, and if the system would relax to a new equi¬ 
librium, the latter would have a magnetization m-eq = 0. 
The stationary state reached after the damping of the os¬ 
cillations is a QSS, very far from the homogeneous equi¬ 
librium at e*. The nonequilibrium character of this state 
is further shown by the fact that the temperatnre profile 


lZodpf{'d,p) 


(5) 


is non-uniform, and there is temperature inversion, as 
shown in Fig. 2, where T{'d) is plotted together with the 
density profile 

/ CO 

dpf{'d,p). (6) 

-OO 


The temperature profile indeed remains essentially the 
same for the whole lifetime of the QSS, as we checked by 
measuring an integrated distance ^ between the actual 
temperature profile and the constant equilibrium one, 
Teq, at the same energy, as follows: 

e(t)= r Im<)-Teq|dz?. (7) 


In Fig. 3, is plotted for systems with different val¬ 
ues of N kicked with the same = 10 at to = 100 for a 


FIG. 3. (Color online) HMF: Time evolution of the distance 
from equilibrium temperature Eq. (7), with N increasing 
from bottom to top: N = 5 x 10^ (red), N = 10® (blue), 
N = 2.5 X 10® (black), N = 5 x 10® (purple). Each curve is 
the average over Ur realizations with Ur ranging from 20 for 
N = 5 X 10® to 10® for iV = 5 X 10®. Inset: ^ as a function of 
t/N. 


duration r = 1. After the kick, ^(t) oscillates and then 
reaches a plateau whose duration grows with iV, as ex¬ 
pected for a QSS. The inset of Fig. 3 shows that if times 
are scaled by N, the curves reach zero at the same time, 
consistently with the lifetime of an inhomogeneous QSS 
being proportional to N [28]. 

We performed simulations starting with equilibrium 
states corresponding to different values of mo and apply¬ 
ing perturbations of different strengths h and duration r, 
and the system almost always ended up in a QSS showing 
at least a partial temperature inversion over some inter¬ 
val of the values of -d, if not the whole; we also considered 
another model where also the second Fourier mode is re¬ 
tained, again obtaining the same qualitative behaviour 
(see Appendix A). 

Preparing the system in equilibrium and then bringing 
it out of equilibrium by means of a perturbation act¬ 
ing for a short time mimics processes that may actually 
happen in nature, e.g., a transient density perturbation 
in a self-gravitating fluid. Hence, as a second example 
we considered a two-dimensional self-gravitating system 
(2DSGS), that is, a system of N particles of mass m 
moving in a plane with Hamiltonian 




N I ,2 

Pil 


i=l 


2m 


N 

- Gm® In 
ip>i 


\/(r,-rj)2-br2 


, ( 8 ) 


where s is a length scale, Pi = {xi,y^), 

tq is a small-scale cutoff and G is the gravitational con¬ 
stant. Such a system can be seen as a simple model of 
filamentary interstellar clouds [29]. We performed MD 
simulations considering = 3 x 10^ particles initially in 
a thermal equilibrium state, whose radial density profile 
is known [30] (see Appendix A), and we kicked them out 
of equilibrium by instantaneously adding to all their ra¬ 
dial velocities the same amount Svr = crr/2, where tr,- is 
the radial velocity dispersion. As in the HMF case, af¬ 
ter the kick the system develops macroscopic oscillations, 
that damp out after a time of order Tdyn = r^\/2/GM 
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FIG. 4. (Color online) 2DSGS: temperature profile T{r) (blue 
squares) and density profile n(r)/no (red circles) measured in 
the QSS at t = 90rdyn while starting from a thermal state 
with uniform temperature To = 0.5, in natural units such 
that r*(t) ~ 7 and rg = 10“®. 


where M is the total mass and r° = = 0) is the (ini¬ 

tial) half-mass radius, and eventually sets in a QSS. The 
latter exhibits temperature inversion up to r > r* (t), as 
shown in Fig. 4. The radial profiles n{r) and T(r) are 
obtained by averaging T(r) and n(r), defined by replac¬ 
ing (i?,p) with (r,p) in Eqs. (5) and (6), over the polar 
angle. 

We have thus shown that temperature inversion arises 
in two different long-range-interacting systems when 
kicked out of equilibrium. Moreover, also if we start 
from a state with a Maxwellian velocity distribution but 
a spatial distribution different from the equilibrium one, 
macroscopic oscillations develop, and after the damping 
of the oscillations, the systems end up in a QSS typi¬ 
cally exhibiting at least a partial temperature inversion 
(see Appendix A). Therefore, this phenomenon does not 
depend on the details of the way the initial state is pre¬ 
pared, provided it is inhomogeneous, and that a collective 
oscillation develops when it is brought out of equilibrium. 

A very general mechanism may then be responsible for 
this phenomenology. In the following, we argue that in¬ 
deed a simple and general explanation can be found. As 
stated above, before collisional effects set in, the dynam¬ 
ics of a long-range interacting system is described by the 
Vlasov equation (2). When a state described by a sta¬ 
tionary solution of the Vlasov equation is perturbed, a 
phenomenon called Landau damping [31] occurs, which 
is a kind of wave-particle interaction responsible for the 
collisionless damping of oscillations induced by the ex¬ 
ternal perturbation. The theory of Landau damping is 
well developed, and, in the case of small perturbations 
of homogeneous states, has been recently put on rigorous 
grounds [32] . The phenomenon occurs also when the per¬ 
turbed state is inhomogeneous [33, 34]. The basic physi¬ 
cal mechanism is the following: The perturbation creates 
a wave in the system, as witnessed by the onset of col¬ 
lective oscillations. Consider for simplicity a monochro¬ 
matic wave. Particles that have a velocity v close to the 
phase velocity Uph of the wave will interact strongly with 
the perturbation, while those with v -C Wph or v ^ ?;ph 
will be essentially unaffected. Particles with velocities 


slightly less than tiph will acquire kinetic energy, while 
those with velocities slightly larger than fph will lose ki¬ 
netic energy. Since for a typical f{v) one has df /dv < 0 
(for positive v, and the opposite for negative v), there 
are more particles that gain energy than those that lose 
energy, so that the wave in effect loses energy, and the 
perturbation is damped. After the interaction with the 
perturbation, a Maxwellian f{v) is modified close to tiph 
(and if the wave is stationary or travels in both directions, 
as in our case, also close to — Uph), where a “shoulder” is 
formed, and (u^) is increased. If, prior to perturbation, 
the state is spatially homogeneous, this happens in the 
same way throughout the system, and the initially coher¬ 
ent energy of the wave goes into uniformly heating the 
system. But in our case the state is not homogeneous 
prior to the perturbation: it is clustered, and the “shoul¬ 
der” in f(v) created by Landau damping is nothing but a 
suprathermal tail, so that now velocity filtration acts and 
produces temperature inversion: fast particles climb the 
potential well higher than slower ones, and suprathermal 
tails grow when density decreases. Indeed, as shown in 
[21], velocity filtration always works in systems described 
by the Vlasov equation (2); there, the suprathermal tails 
of f{v) were fed by the initial conditions, while here it 
is the perturbation that creates them, via Landau damp¬ 
ing. It is interesting to note that standing radiofrequency 
waves are used to heat laboratory plasmas by creating a 
shoulder in the velocity distribution [35]; in that case, 
one needs to sustain the wave from outside, while in our 
example it is velocity filtration that amplifies the effect. 

The above physical picture is idealized and one should 
take into account the coupling of the various modes of the 
perturbation in the inhomogeneous case [33, 34]. This 
notwithstanding, the main point is that the perturbation 
does not interact in the same way with all the particles, 
but preferentially gives energy to particles with already 
rather large velocities. This results in a suprathermal 
velocity distribution, as shown in Fig. 5 (left panel) for 
the HMF model. Soon after the kick a high-p tail shows 
up, which is built up from peaks at different values of 
p corresponding to different positions. The evolution of 
these peaks results in oscillations in /(p) for t < 150, then 
the oscillations damp out and due to the absence of an 
efficient mechanism able to evenly redistribute this excess 
energy among all the particles, the velocity distribution 
stays nonthermal and essentially the same for times t < 
Tcoii, allowing velocity filtration to produce temperature 
inversion. The distribution function of the HMF model 
in the QSS is plotted for two different values of d in Fig. 
5 (right panel), and the growth of the suprathermal tails 
in the less dense parts of the system is well apparent. 
We found similar results for the distribution functions 
also in the 2DSGS case (see Appendix A). For the HMF, 
perturbations of the form cos {kdi) with k ^ 1 

also yield preferential absorption of energy around values 
of V that are different from those in the fc = 1 case, 
coherently with the above picture (see Appendix A). 

We have shown that nonequilibrium stationary states 
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FIG. 5. (Color online) HMF: (left) Momentum distribution f{p) at t — 0 (black dashed line), t = 101 (red plus signs), t = 150 
(magenta empty squares), t = 200 (grey filled triangles), t = 10® (orange empty triangles), and t — 10^ (magenta filled rhombi). 
(Right) Distribution function /(■!?,p) measured in the QSS as in Fig. 2, at = 0 (red squares) and = tt (blue triangles). The 
distribution functions are plotted against p^sgn(p) to better show the difference with respect to the initial Maxwellian. 


with temperature inversions are the typical outcome of a 
perturbation acting for a short time on a clustered equi¬ 
librium state of a long-range interacting system. This 
rather surprising result can be explained in terms of Lan¬ 
dau damping and velocity filtration, suggesting that tem¬ 
perature inversions may occur whenever the dynamics is 
collisionless up to the relevant timescales. This mecha¬ 
nism may be actually relevant to understand temperature 
inversions observed in nature. 
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Appendix A: Supplemental material 

In this Appendix we present some supplemental mate¬ 
rial and in particular, report on the results of molecular 
dynamics simulations performed using different protocols 
and sets of parameters with respect to those shown in the 
paper. Moreover, we show some further details and re¬ 
sults on the two-dimensional self-gravitating system and 
finally we present a tutorial discussion of the mechanism 
of velocity filtration. 


1. HMF model 

Let us start with results for the HMF model, whose 
Hamiltonian is 

N 2 .. N N 

'Hhmf = + (Ai) 

i—1 i—\ j>i 

where di € (—tt, tt] is the angular coordinate of the ith 
particle and pi is the conjugated momentum. The system 
is prepared in a thermal equilibrium state with temper¬ 
ature T and magnetization m = mg, where m is the 


modulus of the vector 

N N \ 

^cosdi,^sindi I . (A2) 

i=l i=l / 

After a time tg >0, an interaction with a magnetic field 
h is switched on by adding to the Hamiltonian a term 

N 

Hh = —h'^cosdi (A3) 

i=l 

and then switched off again at a time tg -}- r. 

a. Density and temperature profiles 

In Fig. 2 of the article, we have shown the density 
profile n{d) and temperature profile T{'d) in the QSS at 
t = 10^, obtained with mg = 0.521 (corresponding to 
T = 0.4244), tg = 100, h = 10, and r = 1. In the 
following, and in particular, in Figs. 6, 7, and 8, we plot 
n and T, still measured at t = 10^ after a kick given at 
tg = 100, but with different values for other parameters, 
as described in the captions. 

In Figs. 6 and 8, we see temperature inversion on the 
whole range of d’s. In Fig. 7, inversion is present only on 
about half the interval, corresponding to —1.5 ^'d ^ 1.5, 
where the state is sufficiently clustered. In the rest of 
the interval, the density is almost flat, and there is no 
temperature inversion, although the temperature profile 
is still nonuniform. 

In Fig. 9, we plot n and T at t = lO'^ as obtained using 
a different protocol. Instead of preparing the system at 
thermal equilibrium, and then adding an impulsive per¬ 
turbation, we prepare the system with a thermal velocity 
distribution at T = 0.4244, but with a spatial distribu¬ 
tion corresponding to a value m of the magnetization that 
it slightly smaller than the equilibrium value mg = 0.521, 
that is, m = 0.516. Again, the system ends up in a QSS 
with temperature inversion, but now the effect is less pro¬ 
nounced, since the state is much closer to an equilibrium 
one than in the other cases. 
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FIG. 6. HMF: Temperature profile r(i?) (blue solid line) and 
density profile n{'&) (red dashed line) measured in the QSS at 
t — 10^, in the same conditions as Fig. 1 of the article, with 
the difference that T — 0.25. 



1 

0.9 

0.8 

0.7 ^ 

0.6 
0.5 
0.4 


FIG. 9. HMF: Temperature profile T{i}) (blue solid line) and 
density profile n(i?) (red dashed line) measured in the QSS at 
t = 10^, while starting from a a Maxwellian velocity distribu¬ 
tion at temperature T = 0.4244, but with a spatial distribu¬ 
tion whose magnetization is smaller than the corresponding 
equilibrium one by an amount Am = 0.05. 


bution into a non-thermal one, due to Landau damping. 
The fact that the energy exchange between the perturba¬ 
tion and the particles occurs via a kind of wave-particle 
interaction is witnessed by the fact that energy is pref¬ 
erentially absorbed by particles having velocities close to 
some given values (see Fig. 4 of the article). As a further 
support to this explanation, in Fig. 10, we show that if 
we change the shape of the perturbation from that given 
in (A3) into 


FIG. 7. HMF: Temperature profile r(i?) (blue solid line) and 
density profile n{'&) (red dashed line) measured in the QSS at 
t = 10^, in the same conditions as Fig. 1 of the article, with 
the difference that h = 1.0. 


b. Momentum distributions 


N 

"H?! = —cos(fcdi) (A4) 

i=l 

with fc ^ 1, we still have preferential absorption of energy 
for specific values of the momentum, but these values are 
different from the case k = 1 and depend on k. 


In the article, it is argued that temperature inversion 
is produced by velocity filtration after the perturbation 
has changed the initially Maxwellian momentum distri- 
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2. Extended HMF model 


In order to show that including further Fourier com¬ 
ponents in the interaction potential does not destroy the 
phenomenon of temperature inversion, we now report on 
some results obtained for an extended HMF model, where 
also the second mode in the Fourier expansion of the po¬ 
tential is retained, defined by the Hamiltonian 


N 2 . N N 

HhMF2 = 

i=l i=l j>i 

- (1 - A:)cos(2i?i - 2'dj)] . 




(A5) 


Note that for AT = 1 one recovers the HMF model. 


FIG. 8. HMF: Temperature profile T{'&) (blue solid line) and 
density profile n(i}) (red dashed line) measured in the QSS at 
t = 10"^, in the same conditions as Fig. 1 of the article, with 
the difference that the duration of the impulse is r = 4. 


a. Density and temperature profiles 

We considered the extended HMF model (A5) with 
K = 0.3. Also in this case, starting from equilibrium 
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FIG. 10. HMF: Momentum distribution at t = 10^, in the 
same conditions as Fig. 1 of the main text, with the difference 
that the impulsive perturbation is of the form (A4) with k = 2 
(top panel) and k = 8 (bottom panel). 
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FIG. 11. Extended HMF model: temperature profile T{'d) 
(blue solid line) and density profile n(r?) (red dashed line) 
measured in the QSS at t = 10®, while starting from a ther¬ 
mal equilibrium state at temperature T — 0.25 kicked out of 
equilibrium by an external held h = 10 acting for one time 
unit. Here N = 10^. 


and kicking the system by switching on an external held 
for a short time we observed temperature inversion. In 
this case the inversion is limited to the region where the 
system is more collapsed, but still the temperature dif¬ 
ference is large (see Fig. 11). 


3. Two-dimensional self-gravitating system 


In the article we have also considered a two- 
dimensional self-gravitating system (2DSGS), that is, a 
system of N point particles of mass m moving in a plane 
with Hamiltonian 


N I 


2m 


N N 

-b Gm^ 

i—1 j>i 


v/(r,-r,)2+r2 


(A6) 


where = {xi,yi), s is a length scale, Pi = (±i,i/i), 
tq is a small-scale cutoff to avoid the divergence when 
\ri — r^l —>■ 0, and G is the gravitational constant. Such 
a system can be seen as a simple model of hlamentary 
interstellar clouds. 


a. Thermal equilibrium state 

The simulations reported in the article have been per¬ 
formed considering N = 3 x 10^ particles initially in a 
thermal equilibrium state, whose radial surface mass den¬ 
sity prohle is [30] 


Q{r) 


go 

(1 + 1 -^)' 


(A7) 


where x = rjs and go is the central surface density. By 
choosing 


s = 


1 


inGPgo 


where 


/? = 


m 


(AS) 


(A9) 


is the parameter fixing the width of the equilibrium 
Maxwell distribution of the velocities and fes is the Boltz¬ 
mann constant, one gets that j3 depends only on the total 
mass (per unit length, since we are considering a two- 
dimensional slice of an infinite cylinder) M of the system 
[30]: 


13 = 


2 

■ 


(AlO) 


Note that M is a finite quantity even if the system is 
unbounded in the plane due to the confining nature of 
the logarithmic two-dimensional gravitational potential. 
We used units where G = M = 1, and we also fixed 
s = 1, so that the equilibrium value of (3 corresponding 
to the prohle (A7) is /3 = 2 and the central density is 
hxed by Eq. (A8), yielding go = (47r)~ and a half-mass 
radius cc* = r*/s = \/8. 
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FIG. 12. 2DSGS: Distribution function f{r,Pr) measured in 
the QSS as in Fig. 3 of the article, at r = 7.9 (red squares) 
and r = 1.4 (blue triangles). 


4. Velocity filtration and temperature inversion 

Velocity filtration was first suggested by Scudder [18- 
20] to explain the heating of the solar corona. We now 
describe the original argument by Scudder in the simple 
case of noninteracting particles in an external potential, 
and then discuss why such a mechanism works also for a 
long-range interacting system, as shown in [21]. 

Let us consider a system of noninteracting particles 
moving in one spatial dimension on a semi-infinite line 
a: > 0. The particles are subjected to an external poten¬ 
tial ip{x), such that the potential energy is growing with 
X, 4’'[x) > 0. A typical case is that of an atmosphere, 
that is, a gas of particles in a gravity field generated by 
a planet. In this case, x is the height above the ground, 
and the potential energy ip^x) is the gravitational one, 


b. Distribution functions 


tfix) 


GM X 
R x + R' 


(All) 


In the article we reported on simulations where, start¬ 
ing from the equilibrium described above, particles were 
kicked out of equilibrium by instantaneously adding to all 
their radial velocities the same amount 5vr = where 
Or is the radial velocity dispersion, that is, Or = 

In those conditions we observed temperature inversion in 
the QSS up to radii larger than r*. In Fig. 12 we plot 
the distribution function f(r,pr) measured in the QSS at 
two radii, that is, r = 1.4 where the temperature is close 
to its minimum and the density is close to its maximum, 
and r = 7.9 where the temperature is close to its max¬ 
imum. The growth of the suprathermal tails in the less 
dense parts of the system as well as “shoulders” around 
|pr| — 1 are well apparent. 


c. Other simulation protocols 


We also observed temperature inversion with other 
simulation protocols. 

We found that temperature inversion also arises in 
long-lasting transients of the dissipationless collapse of 
initially “cold” (i.e., initial virial ratio 2 AT/] V| < 1, where 
K and V are the totale kinetic and potential energies, 
repectively) two-dimensional self-gravitating systems. In 
this case, the initial condition for the numerical simula¬ 
tions was, as before, the density profile given in Eq. (A7), 
while velocities were extracted from a thermal distribu¬ 
tion with P = 10. Radial temperature profiles showing 
inversion up to the radius enclosing 70% of the total mass 
M last for times t ~ 20tdyn, after the violent relaxation 
phase, before the system settles to an isothermal profile 
with an abrupt fall at large radii. Perfectly cold initial 
conditions (i.e., 2A'/| V| = 0 at t = 0) were excluded from 
this study as they are known to be prone to the so-called 
bar instability breaking the system’s axial symmetry (see 
e.g. [4]). 


where G is the gravitational constant, M is the mass of 
the planet, and R its radius. The one-particle distribu¬ 
tion function f(x, p, t) is a solution of the Vlasov equation 


dt ^ dx dx dp 


(A12) 


with a given stationary boundary condition at a; = 0, i.e.. 


Mp) = f{0,P,t) . (A13) 


The Vlasov equation (A12) is nothing but energy con¬ 
servation for each particle, so that due to the potential 
tp{x), there will be a “velocity filtration” effect. This 
means only those particles whose kinetic energy k at 
X = 0 is sufficiently large to overcome the potential bar¬ 
rier Atjj{x') = tp{x') — '0(0) will reach the position x' 
where their kinetic energy will he k' = k — A'0(x'). As a 
consequence, the spatial density 


n{x) 



dpf{x,p) 


(AM) 


is a decreasing function of x. We now consider a 
Maxwellian boundary condition. 


where no = n(0), and for simplicity we have taken parti¬ 
cles of unit mass and set the Boltzmann constant kg to 
unity. Then, the stationary solution f{x,p) of Eq. (A12) 
is the so-called “exponential atmosphere”: 


f{x,p) = exp 


i’jx) - Vi(o) 
To 


/o“(p)- 


(A16) 


The system is isothermal, i.e., the temperature profile 

1 

T{x) = --- / dpp^f{x,p) 
n[x) 


(A17) 
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FIG. 13. Logarithm of the distribution function f{x,p) as 
a function of the signed kinetic energy k calculated for non¬ 
interacting particles in a gravitational potential (All) with 
GM = 10 and i? = 1 at height x = 0 (black dotted line), 
X = 0.25 (red line), and x = 0.65 (blue line). The black dot¬ 
ted line coincides with the boundary condition as given 
by Eq. (A15). 

is constant, T{x) = Tq. In this case, the only effect of 
velocity filtration is that the density decreases with the 
height X. However, this is a very special case. Indeed, 
while velocity filtration occurs for any boundary condi¬ 
tion /o(p), it yields a constant temperature profile only 
if the boundary condition is Maxwellian. This is due to 
the fact that f{x,p) given by Eq. (AI6) is separable in x 
and p. 

In Ref. [18], this was explained using a graphical con¬ 
struction. If we plot liif{x,p) as a function of the kinetic 
energy k = p^ 12, we get a straight line. Now f{x, k) is 
obtained from /(O, k) by removing the part of the dis¬ 
tribution corresponding to kinetic energies smaller than 
the potential barrier A'0(a;), and then rigidly translating 
the remaining part towards the origin. We thus get a 
straight line with the same slope, and since in this case 
T = —(dlnf/dk)~^, the temperature is the same at any 
X. In Fig. 13 we plot In / as a function of the signed 
kinetic energy 

k = sgn(p)y (A18) 

at different heights x in a gravity held given by Eq. (All) 
with GM = 10 and R = 1, calculated solving Eq. (A12) 
with the Maxwellian boundary condition given in Eq. 
(A15). In Fig. 14 we show that the curves collapse onto 
each other if rescaled by the density n{x), as follows from 
Eq. (A16). 

However, the same reasoning implies that if the tails 
of /o(p) are more populated than in a Maxwellian (i.e., 
fo{p) is suprathermal), then velocity filtration yields a 
broader distribution at x than at the boundary. To give 
an example, in Fig. 15, we plot the logarithm of distribu¬ 
tion function divided by the density. In [f{x,p)/n(x)], as 
a function of k at different heights x in the same gravity 
held as in Fig. 13, with a power-law boundary condition 


FIG. 14. Logarithm of the distribution function divided by 
the local density, f{x,p)/n{x), as a function of the signed 
kinetic energy k in the same conditions as Fig. 13. 



k 

FIG. 15. Logarithm of the distribution function normalized 
by the density, f{x,p)/n{x), calculated for noninteracting 
particles in the same gravity field as in Fig. 13, at height 
X = 0 (black dotted line), x = 0.25 (red line), and x = 0.65 
(blue line). The black dotted line coincides with the boundary 
condition fo{p) oc p~^ as given by Eq. (A19). 


The broadening of the distribution as x increases is well 
apparent, at variance with the collapse of the curves in 
the Maxwellian case, implying that the temperature is 
larger at larger heights. Indeed, T grows as n decreases, 
as shown in Fig. 16 where we plot the density and tem¬ 
perature profiles calculated in the same conditions as in 
Fig. 15. 



FIG. 16. Temperature (blue solid line) and density (red 
dashed line) profiles calculated for noninteracting particles 
in a gravity field in the same conditions as Fig. 15. 
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Why do we expect the same phenomenon to also occur 
in a long-range-interacting system in an inhomogeneous 
QSS? Because for times for which the system stays in a 
QSS, its one-particle distribution function f{'d,p,t) still 
obeys a Vlasov equation, 


energy V is 

N N N 

2=1 jai 2=1 

(A22) 

For the HMF model (Al), we have 


df df d{{u) + ^)df 
dt dd dp 


(A20) 


u(i? — d') = 1 — cos(i? — -&'). (A23) 

When present, the external potential term typically de¬ 
scribes the interaction with a magnetic field h, so that it 
is given by 


where we have written d instead of x because we have in 
mind systems with periodic coordinates like the HMF 
model (Al). The only difference with respect to Eq. 
(A12) is that the potential is now composed of two terms: 
the mean-field interaction. 


{u){'d,t) = J M' j dp'u{-d--d')f{d',p',t), (A21) 


where u is the two-body interaction potential, and (pos¬ 
sibly) an external field such that the total potential 


=—hcos'd. (A24) 

Hence, despite the fact that each particle is strongly in¬ 
teracting with all the other particles, the dynamics of 
our system in a QSS is equivalent to that of noninteract¬ 
ing particles immersed in a self-consistent field given by 
(u) + tp. From the above discussion, it should be clear 
that the key ingredients to obtain a temperature inver¬ 
sion due to velocity filtration in the noninteracting parti¬ 
cle case are (i) the presence of an attractive external field, 
and (ii) the fact that /(p) is suprathermal. Hence, we 
should expect the same behavior also in the long-range¬ 
interacting case in a QSS, provided the resulting effective 
field (u) + Ip is attractive—or, equivalently, the density 
profile is clustered due to attractive self-interactions— 
and f{p) is suprathermal. 
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